Revisiting the core dynamics in high-order harmonic generation 

using Bohmian trajectories 
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A full quantum model of high-order harmonic generation is presented from a Bohmian-mechanical 
perspective. According to the three-step model, this phenomenon occurs due to the laser-induced 
recombination of an electron ejected by tunnel ionization with its parent ion. However, when 
revisited within the Bohmian scenario, we find that the high-harmonic spectrum is generated by 
those trajectories that reside well inside the core rather than by those that undergo excursions 
out of it. This agrees with the outcome of the full solution of the time-dependent Schrodinger 
equation, in which the spectrum is obtained through the dipole acceleration. Furthermore, it shows 
that the innermost part of the core, represented by a single Bohmian trajectory, leads to the main 
contributions to the high-harmonic spectra. Nevertheless, one may relate time-frequency maps from 
these central Bohmian trajectories to classical electrons behaving according to the three-step model. 
This happens because the quantum phase carried by each Bohmian trajectory is influenced by the 
whole wavefunction and, therefore, also by those trajectories that leave the core. 

PACS numbers: 32.80.Rm, 42.50.Hz, 42.65.Ky 



High-order harmonic generation (HHG) is a very im- 
portant phenomenon in strong-field physics [J-Q , with a 
wide range of applications, such as subfemtosecond pulses 
[1-0] , and the attosecond imaging of matter (il-Hoj . These 
applications are strongly based on the current paradigm 
for describing HHG: the three-step model (TSM) [lTI,[ia|. 
According to the TSM, HHG takes place when an elec- 
tron, initially ejected from the core by tunneling, comes 
back and recombines with its parent ion, emitting high- 
frequency radiation. This model, in both its classical 
[ill and quantum mechanical formulations, within the 
strong-field approximation (SFA) [12], has become a cor- 
nerstone for theoretical and experimental research on the 
subject. Furthermore, the physical picture of an electron 
returning to and recombining with the core has also been 
inferred from the direct solution of the time-dependent 
Schrodinger equation (TDSE). The HHG time profile, for 
high enough harmonic frequencies, coincides with the re- 
turn times of a classical electron [l3|, [l4| • 

It is also well established that "cleaner" HHG spec- 
tra, with a large plateau and a well-defined cutoff, are 
obtained from the expectation value of the dipole ac- 
celeration a(i) = (*(t)| - W|*(t)), rather than from 
the dipole length in TDSE simulations [H Q]|. The 
dipole acceleration probes regions near the core, while 
the dipole length emphasizes regions closer to the integra- 
tion boundaries [l6j . This suggests that regions near the 
core, where the overlap between continuum and bound 
dynamics is likely to occur, are important to HHG. A 
legitimate questions is, however, whether one can single 
out a specific region in the core as being the most relevant 
to HHG. Apart from that, one may ask the question of 
how the above-stated overlap relates to the physical pic- 
ture propagated by the TSM of an electron undergoing 
long excursions outside the core. 



In this Letter we address these questions using 
Bohmian mechanics [13, EH, where well-defined trajec- 
tories can be assigned to quantum systems without vio- 
lating any of the fundamental postulates and principles 
of standard quantum mechanics. This is possible because 
Bohmian trajectories contain quantum phase informa- 
tion, which classical trajectories do not. Furthermore, 
in contrast to other orbit-based approaches, such as the 
SFA [l2[, the influence of the binding potential is fully 
incorporated. Recently, Bohmian mechanics has been 
applied to strong-field physics at a descriptive and inter- 
pretational leve l |19t — 221 as well as a source for numerical 
algorithms [22|, l23l j. These works essentially follow the 
traditional scheme of considering a set of Bohmian trajec- 
tories and comparing their statistics with the correspond- 
ing quantum results. Here, on the contrary, we focus on 
the information obtained from the analysis of individual 
Bohmian trajectories showing that: (1) the HHG spec- 
trum can be directly obtained from a single Bohmian 
trajectory (no statistics is needed), rather than from the 
dipole acceleration; (2) the specific trajectory that best 
reproduces the HHG spectrum remains in the core, in- 
stead of undergoing excursions out of this region; and 
(3) still, the physical picture related to the TSM can be 
extracted from the time profiles of the above-mentioned 
Bohmian trajectory. This shows that the innermost re- 
gion of the core best reproduces the full HHG spectra. 

Here, we employ a one-dimensional soft-core model, 
which suffices to illustrate the time-dependent electron 
dynamics, and yet contains the essential physical ele- 
ments for linear polarization. Making use of the dipole 
approximation and the length gauge, the electron dynam- 
ics are described by the Hamiltonian (in atomic units) 
H = — | V 2 - ^J, +1 — xEof(t). The time dependence of 
the laser field is given by a flat-top function, 
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FIG. 1: Dipole acceleration (a) computed from the full TDSE solution and its frequency spectrum (b) for a pulse of frequency 
cjo = 0.05 a.u. and field amplitudes Eo — 0.05 a.u. (black solid line) and Eo = 0.085 a.u. (red dotted line). The vertical arrows 
in the lower panels indicate the theoretically predicted positions of the cutoff frequencies, u c « 29.2u;o and ui c ~ 59.2cjo for 
the lower and higher driving-field intensity, respectively. Set of 21 Bohmian trajectories for Eq = 0.05 a.u. (c) and frequency 
spectrum (d) of those with initial positions at x\(0) = —3 a.u. (red dotted line) and xn(0) = a.u. (black solid line). Set of 
21 Bohmian trajectories for Eq = 0.085 a.u. (e) and frequency spectrum (f) of those with initial positions at xg(0) = —0.6 a.u. 
(red dotted line) and xn(0) = a.u. (black solid line). 



(t/r ) sin(cj £), < t < r on 

f(t) = { sm(cj i), r on < t < T off , 

[1 - (i - r off )/r on ] sm(u! t), r off < t < 77 

(1) 

where ujq and To = 2tt/ujo are the frequency and period 
of a field cycle, respectively. r on = 2.25 t is the time 
required for the external field to reach its maximum am- 
plitude, T ff = 2.2570 + ^70 stands for the duration of the 
field with maximum amplitude, and 77 = r Q ff + 2.25 To is 
the time at which the field is switched off (and our quan- 
tum simulations are completed) . We have considered this 
particular type of pulse, with N — 10 (t/ = 14.5t ), to 
ensure a high monochromaticity and therefore an optimal 
resolution of the HHG spectral features; as discussed else- 
where [24], other types of pulses, though more realistic, 
do not provide the same cleanliness for the same param- 
eters used here. We assume that the system is initially 
in the ground state of the unperturbed soft-core poten- 
tial, with energy e = —0.66978 a.u. Following the usual 
scheme in Bohmian mechanics, first the wave function 
is expressed in polar form ^{x,t) — p 1 / 2 {x^t) e %s ^ x ' t '' tl , 
with the probability density p and the phase S being real- 
valued functions of space and time. Then, the Bohmian 
trajectories are obtained [25j after integrating the (also 
real-valued) guidance equation x — VS at each time- 
step, i.e., once *ff(x,t) is known. The latter is obtained 
by exactly solving the TDSE using the split operator 
method. A total of 21 initial conditions have been evenly 
chosen within the interval xq € [— x c ,x c ], with x c = 3 
a.u. In this way, the time evolution of most of the prob- 
ability density will be well monitored, as the integrated 



probability for \x\ > x c is only 0.41% of the total proba- 
bility. 

The dipole acceleration a(t) and its power spectrum, 



a{t) e~ luJt dt 



(2) 



are displayed in Figs.QJa) an( i (b), respectively, for two 
driving-field amplitudes. In both cases, the cutoff fre- 
quencies lie at w c « |eo| + 3.17t/ p , where U p — Eq/4ojq is 
the ponderomotive energy, which agrees with the TSM. 
According to this model, the origin of these two HHG 
spectra is the same: the subsequent recombinations at 
the core (x = 0; x ^ indicates the continuum, i.e., 
an unbound dynamics) of an electron that is ejected via 
tunneling and subsequently re-crosses x = under the ac- 
tion of the external oscillating field. However, when the 
time evolution of the probability density is followed by 
Bohmian trajectories, which play the role of tracer parti- 
cles (2^ |. a very different picture emerges. In particular, 
the concept of tunneling itself, as such, does not exist 
in Bohmian mechanics: similarly to classical argumenta- 
tions, only those trajectories with an appropriate initial 
energy /momentum are able to surmount the correspond- 
ing barrier; otherwise, the trajectories bounce backwards 
27 1. 



In panels (c) and (e) of Fig. [TJ two subsets of 21 tra- 
jectories are displayed illustrating the flow of probabil- 
ity density. For low intensities, they follow back and 
forth the oscillations of the external field. In contrast, at 
higher intensities, while there is a subset of trajectories 
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that displays a similar behavior, part of the probability 
density is being ionized, with the corresponding trajecto- 
ries eventually escaping from the core. In the TSM, the 
bound dynamics are restricted to x — 0. In the present 
model, however, these dynamics essentially seem to span 
the x values between the top of the narrowest effective 
potential barriers V^\x) = —1/^/x 2 + 1 ± xE . These 
boundaries correspond to the classical ionization points 
at maximum field, and occur at x m = ±4.299 a.u. for 
E a = 0.05 a.u., and x m = ±3.198 a.u. for E = 0.085 
a.u. Hence, by observing the probability flow, one may 
view the motion of the Bohmian trajectories oscillating 
around the initial position as "bound" [see panel (c)], 
while those undergoing excursions far away from their 
initial positions and eventually ionize irreversibly may 
be seen as "not bound" [see panel (e)] . This follows from 
the fact that Bohmian trajectories cannot cross through 
the same configuration space point at the same time [28j . 

Based on the above-mentioned criteria, we analyze 
the power spectrum associated with individual Bohmian 
trajectories, i.e., the spectrum obtained when the n-th 
Bohmian trajectory, x n (t), is substituted into Eq. ([2]) in- 
stead of a(t). These spectra have been computed both for 
the lower- and the higher-intensity field [Figs. HJd) and 
(f), respectively]. For the lower intensity, although the 
spectra in panel (d) refer to single trajectories, both con- 
tain information similar to that provided by the dipole 
acceleration [see panel (b)]. Furthermore, the spectra for 
trajectories x\{t) and x\i(t) are very alike, because these 
trajectories display an analogous confined motion. How- 
ever, notice that the closer the initial position of a trajec- 
tory is to x — 0, the more striking the similarity between 
the spectrum obtained for this specific trajectory and the 
dipole acceleration. For the innermost initial conditions, 
apart from the overall intensity, the spectral features of 
a(t) and X\\{t) are nearly identical. This seems rather 
counterintuitive, as (1) one would expect the spectrum 
to arise from the whole set of Bohmian trajectories (prop- 
erly weighted), and (2), according to the TSM, the tra- 
jectories that leave from and return many times to the 
core should play a leading role. 

For the higher-intensity field [Fig. QJf)], the spectra 
of any ionizing trajectory look like that for xg(t), i.e., 
a continuous, monotonically decaying spectrum with no 
harmonic peaks, but only a small signature around the 
fundamental, uj — ujq, and is thus very different from that 
provided by a(t). At first glance, one would be tempted 
to associate these trajectories with the one described by 
the TSM, since the trajectory xg(t) undergoes two rela- 
tively long excursions out of the core before getting even- 
tually ejected. However, once more, it is the innermost 
core trajectory that leads to the correct spectrum, as is 
observed from the analysis of x\\{t). It is worth stress- 
ing that as the initial condition of the Bohmian trajectory 
gets further away from x = 0, not only does the power 
spectrum of the corresponding trajectory lose features of 
the correct HHG spectrum, but it also gains intensity. 
This is precisely the origin of the problems found when 




FIG. 2: (a) Time-evolution of the probability density for 
Eq = 0.05 a.u. and two Bohmian trajectories launched from: 
x = a.u. (BTC) and x = 4.3125 a.u. (BTI). (b) Power spec- 
trum of the two Bohmian trajectories of panel (a). In (c) 
and (d), we display the time-frequency contour-plots for the 
Bohmian trajectories BTC and BTI, respectively, with a win- 
dow function of width a t — 10 a.u. The superimposed lines in 
panel (c) give the classical returning times as a function of the 
electron kinetic energy upon return, for an electron ensemble 
in agreement with the TSM (red dotted line), with xq — and 
Do 7^ (black solid line), and with xq 7^ and vo = (blue 
dashed line). In the latter two cases the soft-core potential 
has been included. For details, see footnote [3lj . 



computing HHG spectra from the dipole length. 

Hence, the HHG spectrum is related to the phase in- 
formation contained in the innermost part of the wave 
function, i.e., that which remains moving inside the core 
for longer times. This phase information, transmitted to 
the Bohmian trajectories through the above-mentioned 
guidance equation, allows us to monitor in an accurate 
fashion the evolution of the wave function [29J and there- 
fore to reproduce the HHG spectrum, just as the dipole 
acceleration does [iEj . 

Now, a legitimate question that arises after these re- 
sults is: If HHG is described by Bohmian trajectories 
close to the core, why is the cutoff energy lu c predicted 
by the TSM recovered? Time-frequency maps constitute 
a widespread tool employed to extract this type of infor- 
mation from TDSE computations [l3l . |30| . These maps 
are obtained through windowed Fourier transforms 

I(w,t) = J s(i)e- ( ' Li),/, 'V' uf 'fl', (3) 

where s(t) denotes the time-series to be analyzed. The 
Gaussian window function introduces time resolution, 
and allows us to blur up interference effects from the 
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spectrum, highlighting classical-like aspects that can be 
compared with the TSM. We have considered the lower 
field intensity case and two Bohmian trajectories: one 
that remains inside the core at any time (BTC) and an- 
other one (BTI) starting in the very-low probability re- 
gion (beyond x c ), which eventually becomes ionized (such 
as in the high- intensity field case). These two trajecto- 
ries, as well as a contour-plot of the probability density, 
are displayed in Fig. [2ja), and their spectra, which fol- 
low the same pattern described above, in Fig.[2Jb). The 
time-frequency maps for these trajectories are shown in 
FigS-HJc) and (d). The map for BTC follows a typically 
periodic structure, in accordance with the classical return 
times predicted by the TSM, while for BTI the map be- 
comes very different [3lj | . Actually, when compared with 
the time-frequency map that one would obtain from the 
dipole acceleration, again we find that in the case of BTC 
the agreement is fairly good [as expected, since the topol- 
ogy of a(t) and xbtc(*) is essentially the same]. 

Summarizing, Bohmian trajectories show that the 
main contribution to the HHG spectrum arises from 
the most internal part of the wave function. It has 
been shown that a single Bohmian trajectory contains 
all the information necessary to obtain the HHG spec- 
trum, namely the trajectory that starts at xq = (and 
x = VS\ t =o = 0), as in the case of the TSM. This is 
a much stronger statement than that provided by the 
dipole acceleration: By using the acceleration, one may 
conclude that the overlap between the continuum and 
bound part of the wavefunction near the core region are 
important. Here, we show that the part of the wavefunc- 
tion located in the immediate vicinity of x = provides 
the HHG spectrum. All other trajectories in the core 
region, when averaged out, contribute to blur up this 
contribution, thus giving rise to a difference in intensity 
between the spectra obtained from the trajectory starting 
at x(0) = and from the dipole acceleration. We have 
chosen several driving-field intensities, frequencies, pulse 



shapes and binding potentials (not only those presented 
in this work) in order to corroborate that this applies 
to regimes with and without significant ionization. At 
first sight, one may have the impression that this con- 
tradicts the TSM. A deeper scrutiny, however, explains 
why the TSM works so well. The TSM reduces the in- 
fluence of the core to a single point, i.e., x — 0. This 
is a good approximation, because the most relevant part 
of ^t(x,t) for HHG is strongly localized. Furthermore, 
quantum mechanically, one cannot attribute "bound" 
and "continuum" dynamics to well defined spatial re- 
gions. Hence, a Bohmian trajectory undergoing confined 
motion may contain both, as shown by the spectral and 
time-frequency analysis of Xn(t). This is possible be- 
cause this Bohmian trajectory evolves under the action 
of the wave function, which not only encompasses local 
information about the space variations of the potential 
function, but also about global changes of the quantum 
phase. This implies that the quantum mechanical phase 
associated with it is influenced by the whole wavefunc- 
tion, i.e., also by the trajectories that ionize and return 
to the core. These latter trajectories, albeit related to a 
very low probability density, affect the quantum phase of 
the core trajectories according to the picture put across 
by the TSM. In contrast, a confined classical trajectory 
does not contain this phase information, and leads to very 
low electron kinetic energies. Hence, classical trajectories 
must revisit the core region recursively in order to probe 
the phase and provide us with the HHG spectrum. 
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